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ABSTRACT 


We  present  a  viscous  computational  fluid  dynamics  (CFD)  simulation  over  two  finite  twisted 
wings  configured  so  as  to  give  a  theoretically  predicted  elliptic  and  parabolic  lift  distributions. 
Local  surface  integration  and  farfield  methods  were  used  to  calculate  the  induced  drag.  The 
objective  of  this  project  is  to  relate  work-potential  losses  (exergy  destruction)  to  the 
aerodynamics  forces  in  an  attempt  to  validate  a  new  design  methodology  based  on  the  second 
law  of  thermodynamics.  Exergy  destruction  for  the  entire  flow  field  was  determined  from 
the  CFD  results.  CFD  results  show  that  the  parabolic  case  produces  smaller  induced  drag 
and  entropy  generation  rates  than  the  elliptic  case.  The  entropy  generation  rates  for  both 
cases  deviated  significantly  from  the  expected  values,  revealing  the  inaccuracy  of  entropy 
generation  rate  prediction  for  a  turbulent  flow.  This  project,  however,  set  up  a  basis  in  terms 
of  analysis  methodology,  from  which  the  future  work  will  follow. 
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FOREWORD 


This  report,  “A  Unified  Methodology  for  Aerospace  Systems  Integration  based  on 
Entropy  and  the  Second  Law  of  Thermodynamics:  Wing  Aerodynamics  Assessment”,  gives  a 
technical  discussion  of  research  accomplished  to  develop  the  methods  to  compute  entropy 
generation  rate  of  different  wing  sections  and  planforms.  It  is  a  part  of  a  larger  effort  to 
develop  the  Second  Law  of  Thermodynamics  into  the  common  currency  for  all  aspects  of 
system  integration.  This  report  considers  the  work  that  needs  to  be  accomplished  by  a  flight 
vehicle  (e.g.,  the  generation  of  lift  to  keep  a  body  aloft  at  sea-level  conditions),  and  identifies 
sources  of  inefficiency  associated  with  finite  wing  aerodynamics. 

The  work  was  led  by  Prof.  Richard  Figliola  of  Clemson  University  guiding  a  team  of 
graduate  students,  Mr.  Jason  Stewart  of  Clemson  University  and  Mr.  Shohei  Nomura  of  the 
University  of  Michigan.  The  Air  Force  Office  of  Scientific  Research  sponsored  the  project 
under  Program  Manager  Dr.  John  Schmisseur.  It  was  administered  by  the  Air  Force 
Research  Laboratory,  Air  Vehicles  Directorate  with  Dr.  Jose  Camberos  as  project  advisor. 

The  work  was  accomplished  during  the  period  24  May  to  20  August  2004. 


Dr.  David  J.  Moorhouse 
Director,  Multidisciplinary  Technology  Center 
AFRL  Air  Vehicles  Directorate 
Wright-Patterson  AFB,  OH  45433 
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SUMMARY 


Using  Air  Force’s  in-house  code,  Air  Vehicle  Unstructured  Flow  Solver,  we 
simulated  the  three  dimensional,  viscous,  laminar  and  turbulent,  low  subsonic  flows  around 
the  finite  wings  with  theoretically  predicted  elliptic  and  parabolic  lift  distributions.  The 
wings  had  a  rectangular  planform  with  geometric  twists.  Using  the  surface  integration  and 
the  Trefftz  plane  methods,  lift  and  the  induced  drag,  respectively,  were  extracted  for  a 
comparison  with  the  classical  lifting-line  theory.  As  expected,  the  results  for  the  lift 
compared  favorably  with  the  theoretical  values;  they  had  errors  of  no  more  than  6%.  The 
induced  drag,  on  the  other  hand,  experienced  a  maximum  error  of  roughly  30%.  We  also 
attempted  to  perform  both  qualitative  and  quantitative  analysis  of  the  entropy  generation 
rates  and  the  exergy  destruction  rates.  For  both  elliptic  and  parabolic  cases,  the  boundary 
layer  seemed  to  produce  the  majority  of  the  entropy.  We  obtained  the  quantitative  figures  by 
taking  the  volume  integral  of  the  total  entropy  generation  rates  over  the  entire  computational 
domain,  including  the  wake.  These  results  deviated  significantly  from  the  approximate 
values,  which  are  simply  the  product  of  the  total  drag  and  the  cruise  velocity. 

Using  the  classical  lifting-line  theory,  we  have  established  the  basic  techniques  for 
building  the  twisted  wing  geometries  with  the  desired  lift  distributions.  The  actual  lift 
distributions  of  these  wings  yet  require  a  validation.  Nonetheless,  we  have  also  established  a 
grid  generation  technique  for  such  twisted  wings.  The  implementation  of  the  farfield  induced 
drag  extraction  method  to  a  viscous  solution  needs  a  further  investigation.  In  particular,  we 
will  examine  the  effect  of  the  location  and  the  size  of  the  Trefftz  plane  on  the  induced  drag. 
As  for  the  entropy  generation  study,  we  will  investigate  the  concept  of  effective  viscosity, 
which  takes  into  account  the  turbulence  effect  lost  in  the  Reynolds  averaging  in  the  Navier- 
Stokes  equations. 

Overall,  the  basis  established  during  the  summer  of  2004  should  provide  a  good 
starting  point  for  the  study  of  an  exergy-based  design  methodology.  The  following  major 
objectives  include  the  grid  refinement  study  and  the  investigation  of  the  relationships 
between  the  aerodynamic  inefficiencies  with  the  irreversible  losses. 


1 


INTRODUCTION 


PROBLEM  DEFINITION 

An  aircraft  represents  an  intricate  system  composed  of  a  number  of  interacting  energy 
and  non-energy  based  sub-systems.  Research  in  aerospace  engineering  has  focused  on  the 
complicated  procedure  of  integrating  such  subsystems  for  some  time.  Historically,  designers 
have  applied  the  design  and  optimization  procedure  only  at  a  subsystem  level.  For  example, 
a  typical  objective  function  for  the  design  of  a  wing  has  sought  to  minimize  drag  for  a  given 
lift  distribution.  Design  of  a  propulsion  system  has  aimed  to  maximize  thrust  while 
maintaining  as  high  thermodynamic  efficiency  as  possible.  Integration  of  these  components 
or  subsystems,  optimized  separately,  requires  a  great  deal  of  experience  and  knowledge.  In 
addition,  this  subsystem-level  optimization  may  ensure  the  optimum  configuration  for  each 
subsystem,  but  does  not  guarantee  an  optimum  design  for  the  entire  vehicle  as  a  whole. 

As  flight  vehicles  evolve,  they  carry  more  and  more  intricate  components  that  require 
distribution  of  power  from  an  energy  source  [1].  With  the  need  to  develop  a  revolutionary 
design  and  optimization  methodology  for  such  systems,  the  basic  approach  started  to  shift 
from  subsystem-level  to  system-level  over  the  last  few  years.  The  main  questions  asked 
include:  Can  we  find  a  single  common  parameter  that  will  hold  across  all  components  in  all 
stages  of  design  process?  To  approach  this  question,  consider  all  the  constituent  components 
as  the  energy-consuming  systems.  One  way  or  the  other,  they  all  have  some  inefficiencies, 
which  can  be  considered  in  terms  of  energy  used  for  work  versus  energy  wasted  [1].  The 
answer  therefore  seems  to  lie  in  this  concept  of  irreversible  losses  that  is  an  additive  quantity. 
So,  can  we  design  a  flight  vehicle  based  not  only  on  the  first  law  of  thermodynamics,  but  also 
on  the  second  law  of  thermodynamics?  We  seek  to  implement  the  second  law  of 
thermodynamics  to  provide  a  common  metric  that  enables  a  comparison  of  the  performance 
of  different  subsystems  in  terms  of  energy  utilization,  and  to  map  in  detail  the  location  of  all 
irreversible  losses.  We  find  a  growing  consensus  in  the  literature  [2, 3, 4, 5, 6]  that  this  should 
lead  to  the  further  advances  of  air  vehicle  design  methodology  and  to  a  true  system-level 
optimization. 

As  part  of  this  system-level  aerospace  vehicle  analysis  and  design  concept,  the  project 
presented  in  this  report  focuses  on  the  aerodynamics  of  low  subsonic  viscous  flow  over  a 
finite  wing.  Flow  over  a  finite  wing  has  been  known  to  produce  number  of  interesting 
phenomena.  Examples  are  tip  vortex  and  trailing  vortices  that  cause  additional  induced  drag. 
A  major  portion  of  this  project  seeks  to  predict  this  induced  drag  and  the  corresponding 
aerodynamic  losses  using  a  3D  viscous  computational  fluid  dynamics  (CFD)  simulation. 

This  task,  however,  does  not  come  so  easy.  In  3D  viscous  flow,  a  variety  of  sources 
subject  the  wing  to  an  assortment  of  aerodynamic  drag.  One  of  these,  skin  friction  drag, 
results  from  the  viscosity  of  air  and  from  velocity  gradients  near  the  wing  surfaces.  Another, 
the  induced  drag  mentioned  above,  results  from  the  generation  of  lift  and  downwash  velocity 
as  generated  by  the  trailing  and  tip  vortices.  Form  drag  constitutes  another  type  of  drag 
resulting  from  an  unbalanced  pressure  distribution  around  the  wing  and  of  viscous  action.  In 
addition  to  all  these  drag  components,  the  inaccuracy  of  a  CFD  calculation  appears  as  a 
spurious  drag.  The  challenging  task  of  extracting  induced  drag  from  a  CFD  calculation 
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comes  from  the  difficulties  encountered  in  numerically  distinguishing  each  drag  source  from 
one  another.  This  research  project  placed  an  emphasis  on  developing  the  induced  drag 
extraction  methods  in  context  of  our  overall  objective  to  apply  the  second  law  of 
thermodynamics . 

Studying  the  energy  efficiency  in  terms  of  second  law  of  thermodynamics  will  help  to 
analyze  the  lost  work-potential  that  can  otherwise  perform  useful  work.  This  thermodynamic 
quantity,  known  as  “exergy”  or  “availability”,  represents  the  maximum  theoretical  work  that 
can  be  obtained  from  a  system  in  taking  it  from  a  given  chemical  composition,  temperature, 
and  pressure  to  a  state  of  chemical,  thermal,  and  mechanical  equilibrium  with  the 
environment  [7],  For  example,  the  fuel  supplied  to  an  aircraft  represents  its  exergy  or 
maximum  work  potential.  To  meet  the  mission  requirements,  exergy  consumption 
(destruction)  must  occur.  Destruction  of  exergy  during  any  chemical  or  physical  process 
equals  the  irreversible  losses  -  entropy  generation  scaled  by  a  reference  temperature.  The 
second  law  of  thermodynamics  suggests  that  we  should  not  necessarily  strive  to  converse 
energy,  but  rather  we  should  attempt  to  conserve  exergy  [8],  In  the  flow  over  a  flight  vehicle, 
minimizing  the  exergy  destruction  saves  the  available  energy  for  distribution  to  other 
components  of  the  aircraft  for  other  types  of  useful  work. 

The  authors  of  this  report  will  study  these  concepts  mentioned  above  using  twisted 
finite  wings  configured  to  produce  two  limiting  (spanwise)  lift  distributions;  elliptic  and 
parabolic.  These  two  special  cases  provide  a  good  starting  point  for  exploring  new  territory, 
since  they  are  known  to  have  special  characteristics.  An  elliptic  lift  distribution  supposedly 
produces  minimum  induced  drag,  according  to  Prandtl’s  classical  lifting-line  theory.  A 
parabolic  lift  distribution  gives  another  special  case  suggested  by  Greene  [9]  to  have 
minimum  entropy  generation.  For  the  elliptic  lift  distribution  in  particular,  classical  theory 
makes  some  major  approximations,  making  it  an  interesting  study.  This  project  could  serve 
as  a  numerical  experiment  to  test  the  limits  of  validity  of  the  well  known  theory.  This  report 
aims  to  present  the  contributions  made  thus  far,  the  current  status,  and  a  number  of 
suggestions  for  future  work. 
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OBJECTIVES 


The  objectives  of  this  project  come  in  several  parts.  Our  first  objective  deals  with  the 
construction  of  twisted  wing  geometries  with  theoretically  predicted  lift  distributions. 
Classical  lifting-line  theory  provides  the  analytical  expression  for  the  required  spanwise  twist 
angle  to  produce  a  desired  lift  distribution.  For  our  second  objective,  we  seek  basic  grid 
generation  techniques  for  constructing  such  twisted  wings.  We  attempt  to  resolve  the  wake, 
as  well  as  the  turbulent  boundary  layer  all  the  way  down  to  the  wing  surface.  These  two 
preparations  lead  to  the  simulation  of  3D,  viscous,  low  subsonic,  incompressible  flow  using 
the  Air  Force  Research  Laboratory’s  in  house  code,  Air  Vehicle  Unstructured  Flow  Solver 
( AVUS)  developed  by  the  Computational  Sciences  Branch  (AFRL/VAAC)  over  the  past 
decade  [11],  Lastly,  the  post-processing  provides  the  visualization  of  the  flow  field,  lift  and 
drag  forces,  and  entropy  generation  and  exergy  destruction  rates.  This  project  should  serve 
as  a  validation  study  for  the  classical  lifting-line  theory  and  a  hypothesis  developed  by 
Greene  [9].  We  also  propose  to  establish  a  solid  basis  for  the  development  of  a  design  and 
optimization  methodology  based  on  the  second  law  of  thermodynamics. 
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COMPUTATIONAL  METHODOLOGY 


We  created  the  wing  geometries  and  the  grid  using  Gridgen  version  15.03  (we  give  a 
brief  description  of  how  the  grids  were  created  below).  The  CFD  code  A  VUS  is  a  cell- 
centered  finite-volume  code.  Its  fundamental  algorithm  utilizes  a  scheme  of  first  order 
accurate  in  space  and  time,  and  the  grid-aligned  exact  Riemann  solution  method  of  Gottlieb 
and  Groth  [11].  The  second  order  accurate  spatial  approximation  method  utilizes  van  Leer’s 
MUSCL  scheme.  For  an  accelerated  convergence,  A  VUS  uses  a  local  time  stepping 
algorithm  for  steady  state  solutions.  However,  the  low  speed  cases  (M  <  0.2)  with  large 
jumps  in  spacing  and/or  cell  size  may  require  time-accurate  calculations  [11],  Because  of 
computational  constraints,  our  grid  suffered  from  large  jumps  in  cell  size,  so,  we  attempted  to 
achieve  a  steady  state  solution  using  the  time  accurate  option.  Implementing  the  time 
accurate  solver  did  not,  however,  yield  fully  converged  solutions.  Poor  grid  quality  in  certain 
regions  may  have  contributed  to  this.  The  computational  constraints  put  a  limit  on  cell  count 
and  consequently  yielded  non-ideal  grids. 

We  attempted  to  obtain  solutions  from  both  laminar  and  turbulent  calculations  with 
the  Wilcox  k-co  model.  Wilcox  k-co  model  is  a  two-equation  Reynolds  Averaged  Navier- 
Stokes  (RANS)  type  model,  which  performs  well  in  the  near-wall  region.  From  these 
simulations,  we  investigated  the  effect  of  turbulence  model  on  the  results. 

The  flight  conditions  used  for  this  study  are  at  a  pressure  of  latm  and  temperature  of 
288K  (sea-level),  at  a  steady  cruise  Mach  number  of  0.2.  The  wing  has  a  rectangular 
planform  with  an  aspect  ratio  of  6  and  a  constant  NACA0012  airfoil  section  of  lm  chord,  and 
a  span  of  6m.  The  wing  has  an  assumed  weight  of  7300  N. 
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GEOMETRY  AND  GRID  GENERATION 


COMPUTATIONAL  DOMAIN 

Since  the  wing  load  and  the  flow  field  are  symmetric  about  mid-span  of  the  wing  for 
a  level  flight,  the  computational  domain  includes  only  a  half  of  the  wing  with  one  boundary 
condition  set  as  the  plane  of  symmetry.  A  total  of  seven  faces  formed  the  computational 
domain  (Figure  4.1).  We  set  all  faces  as  the  far-field  boundary  condition,  except  the  wing 
surface  and  the  face  that  splits  the  mid  span  of  the  wing.  We  labeled  this  mid-span  face  as 
the  symmetry  face  and  modeled  it  with  as  a  solid  wall  with  slip  boundary  conditions.  The 
wing  surfaces  have  no-slip  wall  boundary  conditions.  The  outer  faces  that  enclose  the  wing 
extends  ten  chord  lengths  behind  the  wing,  five  chord  lengths  above  and  below  the  wing,  five 
chord  lengths  in  front  of  the  wing,  and  two  span  lengths  in  width. 


Figure  4. 1 :  Overall  view  of  the  computational  domain 


CLASSICAL  LIFTING-LINE  THEORY 

Our  intent  was  to  construct  a  wing  that  produced  an  elliptic  and  parabolic  lift 
distribution  along  the  spanwise  direction  as  two  separate  cases  to  test  against  classical  lifting¬ 
line  theory.  The  elliptic  lift  distribution  supposedly  produces  minimum  in  induced  drag. 
Greene  [9]  suggests  that  a  parabolic  lift  distribution  produces  minimum  entropy  with  an 
optimal  bending  moment  at  the  root.  Therefore,  we  are  testing  the  hypothesis  that  a  design 
that  minimizes  drag  and  a  design  that  minimizes  entropy  production  are  not  the  same  thing. 
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In  the  classical  lifting-line  theory,  the  lift  distribution  varies  as  a  function  of  the 
effective  angle  of  attack,  the  spanwise  chord  distribution,  and  the  aerodynamic  twist.  In  this 
project,  chord  stays  constant  along  the  span  and  the  wing  has  no  aerodynamic  twist.  Hence, 
the  lift  distribution  becomes  a  function  of  the  effective  angle  of  attack  in  the  spanwise 
direction  only. 

When  a  finite  wing  produces  lift,  its  effective  angle  of  attack  is  less  than  the 
geometric  angle  of  attack,  the  angle  between  the  free  stream  velocity  vector  and  the  chord 
line.  A  high  pressure  region  exists  on  the  lower  surface  relative  to  the  upper  surface  of  the 
wing.  At  the  wing  tip,  the  lift  must  be  zero,  and  therefore  the  upper  and  lower  pressures  must 
equate.  Consequently,  the  flow  at  the  tip  tends  to  move  from  the  high  pressure  region  on  the 
lower  surface  to  the  lower  pressure  region  on  the  top.  This  causes  a  spanwise  flow 
component  along  the  entire  wing  span.  Effectively,  the  wing  suffers  from  a  decreased  lift 
coefficient  and  an  increased  drag  coefficient  relative  to  the  two-dimensional  case  (or  an 
infinite  wing  case). 

This  spanwise  flow  creates  the  wing  tip  vortices  (Figure  4.2  and  4.3).  Since  the  fluid 
on  the  lower  surface  of  the  wing  tends  to  flow  outward  toward  the  tip  and  the  fluid  on  the 
upper  surface  tends  to  flow  inward  toward  the  root,  a  trailing  vortex  forms  when  the  these 
two  flows  meet  at  the  trailing  edge.  The  trailing  vortex  across  the  wing  span  creates  a  vortex 
sheet,  which  induces  a  downwash  velocity  given  by 


w(z)  =  - 


4  n 


bl  2 

I 

-bl  2 


r(z) 


■dz 


(4.1) 


where  T(z)  is  the  circulation  distribution  along  the  span,  Ux  is  the  free  stream  velocity,  and  b 
is  the  total  span  length.  This  downwash  velocity  effectively  reduces  the  approach  angle  of 
the  free  stream  velocity  vector  by  an  amount  of  induced  angle.  The  effective  angle  of  attack 
equals  the  difference  in  the  geometric  and  induced  angles. 


Figure  4.2\  Wing  Tip  Vortex  (Picture  from  DANTEC  Measurement  Technology) 
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Figure  43 :  Wing  tip  vortex  from  CFD  solution 

We  can  achieve  a  particular  lift  distribution  by  giving  the  wing  a  geometric  twist,  an 
aerodynamic  twist,  or  a  combination  of  both.  Geometric  twist  refers  to  a  variation  in  the 
geometric  angle  of  attack  whereas  aerodynamic  twist  refers  to  varying  the  airfoil  sections. 
We  chose  the  geometric  twist  for  this  project  for  the  relative  ease  of  geometry  and  grid 
generation.  The  equation  of  lift  coefficient, 

1  , 

L(z)  =  -m0U„pca0  (z)  (4.2) 


indicates  that  the  effective  angle  of  attack,  «o  must  follow  the  same  trend  as  the  desired  lift 
distribution,  since  the  lift  slope  mo  (2n  according  thin  airfoil  theory),  ffeestream  density  p, 
chord  c,  and  the  freestream  velocity  Ux  stay  constant.  For  instance,  in  order  to  obtain  an 
elliptic  lift  distribution,  ao  must  vary  in  an  elliptic  fashion  as  well.  The  choice  of  the  root 
effective  angle  of  attack  is  mission-specific.  For  an  elliptic  case,  we  chose  an  arbitrary  root 
effective  angle  of  attack  of  5°,  which  corresponds  to  an  assumed  aircraft  weight  of  about 
7300N. 


The  induced  angle  refers  to  the  change  in  approach  angle  experienced  by  the  free 
stream  velocity  due  to  the  downwash  velocity  discussed  earlier.  From  geometry  (Figure  4.4), 
downwash  angle  (a,)  is  simply 


ai(z)  =  tan 


(4.3) 


which  simplifies  to 


«<(z)  = 


M.z) 

u„ 


(4.4) 
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Free  Stream 


for  a  small  downwash  velocity  compared  to  the  freestream  velocity.  For  an  elliptic  lift 
distribution,  the  induced  angle  of  attack  stays  constant  along  the  span: 


a<  = 


2bU 


(4.5) 


Here,  r0is  the  root  circulation,  determined  from  the  5°  root  effective  angle  by 

_  mycUjXy 
0  2 

Now,  the  geometric  angle  as  a  function  of  span  can  be  determined  by 

a(z)  =  a0(z)  +  ai(z ) 


(4.6) 


(4.7) 


The  resulting  function  for  the  geometric  angle  of  for  an  elliptic  lift  distribution  with  a  5° 
effective  root  angle  is  then 


<*(Z)  =  awat 


(4.8) 


where  aroot  is  a  root  geometric  angle  of  attack  of  6.30899°. 

In  order  to  have  an  accurate  comparison,  both  elliptic  and  parabolic  wings  must 
produce  the  same  total  amount  of  lift.  This  requires  that  the  area  under  a  curve  of  effective 
angle  of  attack  vs.  wing  span  should  be  equal  when  both  wings  are  flying  under  the  same 
flow  conditions.  Therefore, 


Solving  this  equation  yields  the  root  effective  angle  of  5.89049°  for  the  parabolic  case.  The 
induced  angle  of  attack  as  a  function  of  span  for  the  parabolic  case  turns  out  to  be 
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(4.10) 


ai(z)  =  - 


2F0 

b2KUx 
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This  equation  indicates  that  the  induced  angle  of  attack  approaches  negative  infinity  as  z 
approaches  bl 2  (wing  tip).  Obviously,  this  presents  a  problem  when  building  the  wing 
geometry.  To  overcome  this  problem  the  slope  of  the  geometric  angle  of  attack  curve  at 
99.7  %  of  the  span  was  used  to  extrapolate  out  to  the  tip.  This  resulted  in  a  root  geometric 
angle  of  attack  of  7.85398°  and  a  tip  geometric  angle  of  attack  of  -4.99262°  for  an  overall 
geometric  twist  of  12.85°.  Interestingly,  the  parabolic  case  has  the  negative  tip  geometric 
angle  of  attack.  Nonetheless,  both  effective  angles  at  the  wing  tip  should  be  0°. 


CONSTRUCTING  THE  WING  GEOMETRY 

The  wings  were  constructed  using  a  straight  leading  edge  that  has  no  slope  in  any 
direction.  We  attached  the  airfoil  sections  to  the  leading  edge,  rotated  the  proper  amount  to 
give  the  desired  geometric  angle  of  attack.  Once  we  placed  an  appropriate  number  of  airfoils, 
an  interpolating  scheme  in  Gridgen  was  used  to  create  the  surface  between  the  airfoil 
sections.  Figure  4.5  shows  an  example  of  the  wing  geometry. 


Figure  4.5 :  Twisted  Wing  Geometry 


GRID  GENERATION 

Constructing  a  proper  grid  perhaps  takes  up  the  most  essential  part  of  a  CFD 
simulation.  Grid  construction  unfortunately  remains  more  of  an  art  than  a  science  and 
consequently  requires  a  significant  amount  of  time  to  develop  a  proper  grid  to  adequately 
resolve  the  essential  flow  physics.  Therefore,  studying  the  flow  field  a  priori  is  essential, 
especially  when  the  solution-adaptive  grid  refinement  is  not  available.  Overall,  the  grid  is 
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hybrid,  consisting  of  the  structured  hexahedral  cells  in  the  boundary  layer,  and  the 
unstructured  prism-shaped  cells  elsewhere. 

In  a  viscous  solution  one  obvious  place  that  requires  a  fine  grid  resolution  is  in  the 
boundary  layers.  There  are  some  guidelines  to  ensure  the  adequate  resolution  in  this  region. 
For  the  flow  conditions  for  this  study,  the  boundary  layer  will  experience  a  transition  from 
laminar  to  turbulent  at  roughly  ten  percent  of  the  chord  from  the  leading  edge.  Because  of 
the  very  high  gradients  at  the  surface  of  the  wing,  especially  for  turbulent  regions,  it  is 
necessary  to  cluster  the  cells  in  that  region. 

In  a  turbulent  boundary  layer  there  exists  a  laminar  viscous  sublayer  that  extends 
from  the  wall  out  to  y+  (a  normalized  distance  that  is  a  function  of  shear  stress)  of  about  five. 
The  first  node  distance  from  the  wall  should  have  a  y+  value  of  one  or  less  when  using  a  two- 
equation  turbulence  model.  This  ensures  that  at  least  a  few  cells  will  describe  the  viscous 
sublayer.  Determination  of  the  first  node  distance  from  the  wall  starts  from  the  flat  plate 
boundary  layer  approximation.  First  the  Reynolds  number  based  on  chord  length  must  be 
determined: 

Re  =  —2—  (4.11) 

o 


This  Reynolds  number  automatically  determines  the  skin  friction  coefficient  with  the 
empirical  correlation  (for  turbulent  flow): 


0.0592 


cf=- 


Re 


0.2 


(4.412) 


Using  this  value  for  the  skin  friction  coefficient,  the  shear  stress  at  the  wall  is  then 


T 


W 


(4.13) 


With  this  shear  stress  approximation,  the  friction  velocity  can  be  found: 


* 

u 


(4.14) 


Now  the  first  node  distance  is  obtained  by  setting  y+  equal  to  one  in  the  equation  relating 
friction  velocity  to  y+: 

u/_  (4.15) 

y  * 

u 

Using  above  relations,  we  see  that  the  necessary  first  node  distance  is  roughly  5x1  O'6 
m.  A  problem  arises  with  such  a  small  length  is  used  to  construct  a  cell  since  the  volume  is 
now  about  10"18m.  This  may  increase  truncation  error  significantly,  and  the  computer  may 
easily  treat  this  volume  as  practically  zero.  In  fact,  the  CFD  utility  code  used  to  check  grid 
quality  for  A  VUS,  Blacksmith,  flags  such  a  grid  as  unacceptable  if  this  first  node  distance  is 
used.  To  overcome  this  problem,  we  scaled  up  the  entire  model  100  times  and  used  the  CGS 
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(centimeter,  gram,  second)  units  when  executing  AVUS.  Figure  4.6  shows  the  resulting  y+ 
contour  plots  on  the  wing  surface. 


Figure  4. 6:  Normalized  distance,  r  contours  on  the  surface  of  the  wing. 


It  has  been  suggested  that  roughly  150  to  200  points  chordwise  on  both  the  top  and 
bottom  of  the  wing  should  suffice  for  an  accurate  drag  prediction  when  used  in  conjunction 
with  the  other  parameters  concerning  the  boundary  layer  grid  [10].  The  leading  and  trailing 
edge  require  clustering  of  the  grids.  This  allows  the  grid  to  accurately  resemble  the  smooth 
curvature  and  provides  a  good  transition  between  cells  on  the  trailing  edge  and  those  in  the 
near-wake.  The  first  node  distance  of  0.450mm  on  the  leading  edge  was  taken  from  [10]. 
We  also  rounded  the  sharp  trailing  edge  to  facilitate  the  grid  generation  on  the  wing  tip  and 
the  near-wake.  To  accurately  trace  this  small  curvature,  a  first  node  distance  of  0.006  mm 
was  used  on  the  trailing  edge. 

Within  the  high-gradient  regions  like  the  boundary  layer,  the  growth  rate  of  the  cells 
also  plays  an  important  role.  It  is  suggested  that  the  growth  rate  should  at  most  be  20-30% 
[10]  to  avoid  a  significant  error  associated  with  the  gradient  calculations.  Another  generally 
accepted  rule-of-thumb  is  to  have  at  least  ten  cells  within  the  boundary  layer  normal  to  the 
wall.  We  approximated  the  boundary  layer  thickness  at  a  tenth  of  the  chord,  using  the 
Blasius  solution  for  the  flat  plate  laminar  boundary  layer: 


8  =  5 


(4.16) 


The  flat  plate  laminar  boundary  layer  solution  was  used  here  since  this  is  approximately  the 
transition  location  from  laminar  to  turbulence.  Using  this  analytical  solution  gives  a  normal 
distance  in  which  10  cells  should  be  placed.  However,  when  considering  the  first  node 
distance  and  a  growth  rate  of  1.2,  there  will  be  more  than  10  cells  in  the  boundary  layer  at 
this  location.  Even  with  the  turbulence  models  assuming  a  turbulent  boundary  layer  over  the 
entire  wing,  the  grid-stretching  scheme  ensures  10  cells  within  the  boundary  layer. 

Another  important  parameter  in  the  boundary  layer  grid  generation  is  the  number  of 
nodes  along  the  span  of  the  wing.  It  was  suggested  through  a  personal  correspondence  [10] 
that  100  points  spanwise  would  be  a  reasonable  number  to  start  with.  However,  since  the 
grid  points  need  to  be  clustered  at  the  wing  tip,  the  number  of  required  grid  points  could  rise 
to  about  200.  The  reason  for  the  clustering  at  the  wing  tip  is  to  allow  for  a  good  spanwise 
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transition  from  a  wing  tip  boundary  layer  grid.  In  this  project,  however,  we  could  only  afford 
eighty  nodes  due  to  computer  resource  limitations.  A  view  of  the  surface  grid  on  the  wing  is 
shown  in  Figure  4.7  and  4.8. 


Figure  4. 7:  Wing  Surface  Grid 


Tight  clustering 
normal  to  wing 
in  boundary 


Figure  4.8:  Wing  Surface  Grid  and  Surrounding 


The  wing  tip  presents  a  difficulty  when  producing  a  high  quality  structured  grid. 
Having  a  sharp  trailing  edge  results  in  highly  skewed  cells  that  could  adversely  affect 
convergence  of  the  solution.  Therefore,  we  rounded  off  the  trailing  edge  at  99.75%  of  the 
chord  to  facilitate  the  grid  generation  procedure.  Since  the  wing  tip  region  will  contain 
complicated  flow,  the  boundary  layer  mesh  and  the  wing  tip  mesh  should  have  as  smooth 
transition  in  size  as  possible.  An  illustration  of  this  is  shown  in  Figure  4.9.  Figure  4.10  and 
Figure  4.11  show  the  overall  wing  tip  grid  and  trailing  edge  tip  grid,  respectively.  The  wing 
tip  was  divided  into  several  sections  to  create  the  grid  shown  in  Figure  4.10.  Again,  an 
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affordable  grid  count  was  limited  by  computer  resources.  Therefore,  the  grid  suffered  from 
discontinuous  jumps  in  size,  which  are  apparent  in  Figure  4.1 1. 


Figure  4.10:  Overall  Wing  Tip  Mesh 


Figure  4. 11\  Trailing  edge  tip  grid 
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The  trailing  vortices  in  the  wake  can  be  resolved  effectively  with  the  unstructured 
grid  without  increasing  the  cell  count  significantly  (Figure  4.12).  Here,  the  grid  is  clustered 
in  the  wake  and  grows  quickly  to  the  outer  region.  We  approximated  the  growth  rate  of  the 
wake  height  as  2x1/3  [13]  where  x  is  the  streamwise  location. 


Figure  4.12\  Slice  through  domain  showing  grid 

A  structured  cylindrical  grid  was  also  constructed  around  the  wing  on  top  of  the  structured 
boundary  layer  grid.  The  surface  of  the  boundary  layer  grid  follows  the  twist  of  the  wing  in 
order  to  maintain  a  constant  node  distance  from  the  surface  of  the  wing.  This  does  not 
provide  a  uniform  surface  to  sweep  a  face  along  to  create  a  volume  grid  for  the  outer  region. 
Therefore,  we  placed  a  cylindrically  shaped  structured  grid  to  provide  a  uniform  surface.  A 
close  up  view  of  this  grid  is  shown  in  Figure  4.13. 


Figure  4.13 :  Cylindrical  grid  around  the  wing 


Once  we  generated  the  grid  for  the  half  of  the  domain  containing  the  wing,  the 
volume  grid  was  completed  by  simply  sweeping  all  faces  that  lie  in  the  mid  plane  of  the 
domain,  including  the  wing  tip  grid,  in  a  linear  fashion. 
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RESULTS 


Following  accepted  guidelines  for  the  grid  generation  resulted  in  a  total  cell  count  in 
the  order  of  several  million.  However,  slow  computational  speed  hindered  the  calculations 
from  reaching  full  convergence  within  the  given  time  limit.  The  machines  used  for  all  the 
calculations  consisted  of  a  Linux-based  Beowulf  cluster  of  Pentium  III  800MHz  processors 
each  with  1GB  RAM.  The  cluster  contained  sixty  four  processors  in  total  and  we  used  the 
maximum  of  thirty  processors  for  any  given  job.  Executing  A  VUS  on  thirty  processors  to 
solve  incompressible  viscous  flow  on  a  grid  with  3.6  million  cells,  for  example,  took  roughly 
thirty  minutes  per  iteration.  Reasonable  convergence  generally  required  at  least  5000 
iterations.  Therefore,  to  achieve  satisfactory  convergence  at  this  rate  would  require  nearly 
104  days. 

To  obtain  results  in  a  reasonable  time,  it  was  necessary  to  coarsen  the  grid.  In 
addition,  we  increased  the  growth  rate  of  the  cell  size  as  much  as  possible  and  coarsened  the 
wake  region  as  well.  We  filled  the  outer  region  with  a  very  coarse  grid  to  reduce  the  number 
of  cells  as  much  as  possible,  so  the  growth  rate  in  this  region  exceeds  1.3.  After  this  grid 
coarsening  process,  the  total  cell  count  dropped  down  to  1.79  and  1.87  million  for  the  elliptic 
and  parabolic  wing,  respectively.  With  cell  count  of  1.87  million,  for  example,  it  took 
roughly  four  days  to  complete  3800  iterations. 

For  every  simulation,  we  let  twenty  processors  work  in  parallel  to  solve  the  flow  as 
close  to  full  convergence  as  possible  to  a  steady  state  under  low  subsonic  condition  (M=  0.2, 
7?ec=4,400,000  and  standard  temperature  and  pressure).  Residuals  for  continuity,  y+,  and 
forces  in  the  x  and  y  direction  respectively,  Fx  and  Fy,  were  monitored  during  the  calculation 
to  check  the  convergence  status.  However,  grid  quality  hampered  the  calculations  from 
reaching  full  convergence.  Therefore,  we  attempted  to  bring  the  solutions  as  close  as 
possible  to  full  convergence.  On  average,  it  took  roughly  5000  iterations  for  y+  residuals  to 
level  off,  but  the  force  components  in  the  x-direction  continued  to  exhibit  large  oscillations, 
even  after  10,000  iterations. 


VALIDATION  OF  GRID  QUALITY 

To  verify  grid  resolution  one  may  check  the  y+  values  on  the  wing  surfaces,  which 
should  be  less  than  unity.  Contours  of  y+  showed  that  for  both  elliptic  and  parabolic  cases,  y+ 
values  stayed  well  below  one  on  all  wing  surfaces,  except  at  the  very  small  region  towards 
the  tip  (Figure  6.1),  where  it  exceeded  two.  Although  these  high  y+  values  lie  in  such  a  tiny 
strip  of  surface  that  they  may  not  be  affecting  the  lift  and  drag  values  significantly,  we 
suggest  that  the  first  grid  size  normal  to  the  wing  surface  in  these  areas  should  be  reduced. 

Figure  6.2  shows  the  pressure  contour  around  the  sectional  plane  at  z  =  2.9m  for  the 
elliptic  case.  Non-smooth  contour  lines  in  the  outer  regions  indicate  the  under-resolved  flow. 
These  regions  require  further  grid  refinement  in  future  work,  since  they  could  affect  the 
solution  in  the  vicinity  of  the  wing  surface. 
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Figure  6.1 :  Region  of  high  y  values 


COMPARISON  OF  LIFT  AND  DRAG  WITH  LIFTING-LINE  THEORY 

One  of  the  objectives  of  this  project  is  to  predict  induced  drag  and  compare  with  the 
theoretical  values  from  lifting-line  theory.  However,  prediction  of  aerodynamic  forces, 
especially  drag  prediction  is  a  very  difficult  task.  Therefore,  we  are  in  a  sense  testing  the 
capabilities  of  aerodynamic  force  prediction  of  CFD. 

There  exist  several  different  methods  for  predicting  aerodynamic  forces.  A  typical 
method  of  calculating  total  drag  and  lift  forces  is  the  surface  integration  technique,  in  which 


17 


the  pressure  and  the  skin  friction  are  integrated  over  the  solid  surfaces.  Often  times,  using 
this  method  yields  an  accurate  lift  prediction,  but  faces  difficulties  in  total  drag  prediction. 
This  inaccuracy  comes  from  the  need  to  approximate  the  smooth  curvature  of  the  object  with 
flat  cells  and  to  numerically  approximate  high  gradients  [14]  in  the  boundary  layer.  These 
approximations  lead  to  production  of  spurious  contributions  to  drag  as  well,  which  appears  in 
pressure  and  skin  friction  over  a  body  surface.  Therefore,  spurious  contribution  cannot  be 
differentiated  from  the  physical  drag  by  the  surface  integration  method  [15]. 

We  therefore  implemented  the  farfield  analysis  of  Trefftz  plane  for  the  induced  drag 
extraction.  Formulating  momentum  conservation  equations  for  the  control  volume  that 
encloses  the  entire  wing  leads  to  the  following  simple  integral  to  be  taken  over  the  cross-flow 
plane  (Trefftz  plane)  at  a  downstream  location  of  the  wing  (Figure  6.3): 

A  =  —P^.fv2  +  w2dS  (6.1) 

Here,  v  and  w  are  the  lateral  velocity  perturbations.  This  integral  simply  sums  the  lateral 
kinetic  energy  across  the  Trefftz  plane.  The  idea  is  that  this  lateral  kinetic  energy  is  caused 
purely  by  trailing  vortices,  the  reversible  phenomena  that  are  the  source  of  induced  drag. 
This  equation  assumes  the  negligible  streamwise  velocity  perturbations  at  the  location  of  the 
Trefftz  plane. 


Figure  6.3:  Schematic  of  Trefftz  plane 


Many  researchers  have  applied  Trefftz  plane  method  to  the  inviscid  case  frequently 
and  obtained  fairly  accurate  drag  and  lift  results.  In  an  inviscid  case,  the  location  of  the 
Trefftz  plane  is  not  important,  as  long  as  the  streamwise  velocity  perturbation  can  be 
considered  negligible,  and  the  plane  encompasses  the  entire  laterally-perturbed  area. 
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However,  when  it  is  applied  to  the  viscous  case,  computation  gets  more  difficult  because  of 
the  natural  viscosity  of  the  fluid  that  tends  to  damp  out  these  lateral  motions.  In  addition,  the 
artificial  viscosity  embedded  in  the  numerical  scheme  itself  also  contributes  to  the  dissipation 
of  the  solution  in  the  wake,  especially  where  the  grid  gets  coarse.  Therefore,  the  location  of 
Trefftz  plane  should  be  where  the  flow  experiences  as  little  influence  from  the  artificial 
viscosity  as  possible.  In  this  project,  we  placed  it  at  six  chord-lengths  downstream  of  the 
wing. 

Table  6.1  shows  the  result  for  the  induced  drag  coefficient  from  the  lifting-line  theory, 
as  well  as  the  results  from  the  CFD  calculation.  Total  drag  coefficients  from  the  surface 
integration  method  are  also  tabulated  for  comparison.  It  is  apparent  that  the  CFD  calculation 
under-predicted  the  theoretical  values.  This  may  be  due  to  the  dissipating  transverse  kinetic 
energy,  which  is  a  consequence  of  both  natural  and  artificial  viscosity.  It  also  predicted  the 
induced  drag  for  the  parabolic  cases  to  be  smaller  than  the  elliptic  cases.  This  contradicts 
with  the  classical  theory  that  says  elliptic  lift  distribution  produces  minimum  induced  drag. 

Table  6.2  shows  the  result  for  lift  coefficient  calculated  by  surface  integration  over 
the  wing  surface.  Elliptic  and  parabolic  cases  have  the  same  theoretical  values,  since  the 
wings  were  designed  to  have  the  same  total  lift  and  they  have  the  same  wing  surface  area. 
Lift  coefficient  values  are  more  accurate  than  the  induced  drag  values. 

One  important  point  about  the  classical  lifting-line  theory  is  that  the  derivation  of 
minimum  induced  drag  for  elliptic  lift  distribution  assumes  a  planar  vortex  sheet  that  extends 
to  infinity  (inviscid  approximation).  In  actual  flow  over  a  finite  wing  however,  the  vortex 
sheet  rolls  up  out-of-plane  (Figure  6.4  and  6.5)  and  dissipates  because  of  viscous  effect. 
Although  classical  lifting-line  theory  is  shown  to  provide  surprisingly  accurate  induced  drag 
distribution  [9],  we  should  not  exclude  this  approximation  from  the  list  of  possible  reasons 
for  differences  when  comparing  theoretical  and  predicted  induced  drag. 

ENTROPY  GENERATION  AND  EXERGY  DESTRUCTION  RATE 

The  most  important  objective  of  this  project  was  to  test  the  capabilities  of  CFD  to 
calculate  the  entropy  generation  and  the  exergy  destruction  rates  and  relate  them  to  the 
aerodynamic  terms,  such  as  lift  and  drag.  AVUS  calculates  the  total  entropy  generation  rate 
in  every  cell  for  viscous  flow  by 
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(6.2) 


•  =  t9  C)ui  q,  QT 
8en  T  dxj  T2  dxj 

and  we  normalized  it  as  follows: 


(6.3) 

In  equation  6.2,  the  first  and  the  second  term  corresponds  to  the  viscous  and  heat  entropy 
contribution  to  the  total  entropy  generation  rate.  For  low  subsonic  incompressible  flow, 
temperature  gradients  will  be  minimal.  Therefore  the  majority  of  the  contributions  will  come 
from  the  viscous  term.  Temperature  gradients  are  calculated  indirectly  from  density  and 
pressure. 

Figure  6.5  -  6.10  are  the  contour  plots  of  the  total  entropy  generation  rate  for  the 
elliptic  turbulent  case.  As  expected,  the  boundary  layer,  leading  edge  region  and  the  near 
wake  produces  the  highest  entropy  generation  rate  (Figure  6.8  and  6.9).  Entropy  generation 
rate  in  the  region  indicated  by  blue  is  essentially  zero. 

Figure  6.11  -  6.14  show  the  contour  plots  of  total  entropy  generation  rate  for  the 
parabolic  turbulent  case.  This  wing  produced  a  similar  total  entropy  generation  rate  plot  as 
compared  to  the  elliptic  case. 

Although  majority  of  the  entropy  generation  happens  within  the  boundary  layer 
where  the  high  velocity  gradients  exist,  it  is  natural  to  think  that  the  wake  and  its  entropy 
generation  is  a  consequence  of  the  moving  wing  disturbing  the  flow.  It  is  also  natural  to 
think  that  wings  with  different  lift  distributions  produce  wakes  with  different  levels  of 
entropy  generation  and  there  may  be  a  connection  between  entropy  generation  rates  over  the 
wing  surface  and  in  the  wake.  Therefore,  to  account  for  both  cause  and  effect,  volume 
integrals  of  the  total  entropy  generation  rates  were  computed  during  post-processing: 

gen) total  ~  \\\s ger>dV  (6.4) 

V 

It  is  interesting  to  note  that  for  steady  motion  the  theoretical  value  of  this  quantity  should  be 
related  to  the  power  consumed  by  the  moving  wing  and  the  reference  temperature  by  the 
following  simple  expression: 
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(6.5) 


(S  ) 

V  sen  )totai  T 

J-  n 


This  equation  simply  says  that  the  power  consumed  by  the  horizontally  moving, 
unaccelerated  wing,  DU.r  should  equal  to  the  exergy  destruction,  T0  (Sgeii )  .  Total  drag,  D, 

was  approximated  from  the  total  drag  coefficient  directly  from  AVUS.  Therefore  equation 
6.5  only  provides  an  approximate  value  to  which  the  results  can  be  compared.  Table  6.4 
tabulates  both  the  reference  values  and  the  CFD  result. 

Notice  the  significant  differences  between  the  expected  values  and  the  results  from 
CFD.  This  could  be  due  to  several  reasons:  First  of  all,  according  to  the  flat  plate  boundary 
layer  theory,  this  problem  should  contain  laminar  region  up  to  roughly  11%  of  the  chord,  and 
transition  and  turbulent  region  elsewhere.  Solving  such  flow  with  laminar  Navier-Stokes 
equations  assumes  laminar  flow  everywhere  in  the  computational  domain  and  therefore  loses 
the  effect  of  turbulence.  Secondly,  the  equation  for  total  entropy  generation  rate  (Eq.6.2)  only 
accounts  for  the  molecular  viscosity.  In  CFD,  Reynolds  Averaging  Navier-Stokes  (RANS) 
types  of  equations  tend  to  dissipate  a  lot  of  small  scale  structures  in  the  flow.  This  leads  to 
simplified  flow  field  and  gradients,  and  apparently  inaccurate  entropy  generation  predictions. 
There  exists  a  concept  of  effective  viscosity  in  turbulence  models  for  RANS  solvers.  Its 
purpose  is  to  recover  the  lost  effect  of  those  small-scale  structures.  The  effective  viscosity  is 
known  to  be  up  to  two  orders  of  magnitude  greater  than  the  absolute  viscosity.  Since  the 
effective  viscosity  value  was  not  implemented  in  these  entropy  computations  (equation  6.2), 
this  important  concept  will  be  investigated  in  the  future  work. 
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DISCUSSION 


As  noted  earlier,  computational  and  time  constraints  limited  our  study  significantly. 
This  led  to  under-refined  areas  within  the  grid  and  undesirable  growth  rates  in  the  outer 
region  as  well  as  very  large  cells  near  the  leading  edge  of  the  wing.  These  factors  should  be 
kept  in  mind  when  examining  the  results  of  this  study. 

As  mentioned  in  the  results  section,  the  lift  values  from  the  CFD  simulations 
compared  favorably  with  the  theoretical  values  from  the  lifting  line  theory.  The  induced  drag 
values  were  calculated  to  be  lower  in  both  cases  as  compared  to  the  theoretical  values.  The 
parabolic  case  resulted  in  the  lower  induced  drag  than  the  elliptic  case,  which  disagrees  with 
the  classical  lifting  line  theory. 

Induced  drag  is  not  very  easily  extracted  from  a  CFD  simulation.  Surface  integration 
only  yields  the  total  drag  from  the  pressure  and  the  skin  friction.  Thus  it  is  not  possible  to 
extract  only  the  induced  drag  contribution  from  the  surface  integration.  Trefftz  plane  method, 
on  the  other  hand,  attempts  to  calculate  the  induced  drag  from  the  lateral  velocity 
perturbations  caused  only  by  the  source  of  induced  drag.  The  difficulties,  however,  of 
applying  this  method  to  a  viscous  CFD  results  were  demonstrated  in  this  project.  As 
expected,  the  calculated  induced  drag  was  smaller  than  the  theoretical  values,  possibly  due  to 
dissipated  lateral  velocity  perturbations.  We  expected  that  the  induced  drag  will  further 
decrease  with  the  distance  from  the  wing.  This  is  indeed  the  case  as  seen  in  Figure  7.1, 
which  shows  how  the  induced  drag  coefficient  decreases  as  the  Trefftz  plane  is  moved 
further  downstream.  Figure  7.1  also  shows  the  effect  that  the  area  of  the  Trefftz  plane  has  on 
the  induced  drag  calculations.  All  Trefftz  planes  extended  six  meters  in  the  spanwise 
direction  but  were  varied  in  the  vertical  direction.  This  was  done  to  see  the  effect  the  area  of 
the  Trefftz  plane  has  on  the  calculations.  The  Trefftz  plane  should  not  extend  to  the  very 
course  outer  regions  to  avoid  the  contamination  from  the  large  cells  and  the  boundary 
condition.  Therefore  an  appropriate  sized  Trefftz  plane  needs  to  be  found.  It  is  apparent 
from  the  results  shown  in  Figure  7. 1  that  a  future  study  needs  to  investigate  this  further  to 
find  a  good  balance  in  the  Trefftz  plane  area. 
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The  steep  drop  in  the  induced  drag  coefficients  at  small  streamwise  distances  in 
Figure  7.1  can  be  attributed  to  the  fact  that  the  streamwise  velocity  perturbation  is  not 
relatively  small  close  to  the  wing.  This  negates  the  use  of  the  integral  of  the  lateral  velocity 
perturbations,  which  assumes  the  negligible  streamwise  velocity  perturbation.  A  modified 
Trefftz  plane  approach  should  be  used  in  these  regions  in  which  the  streamwise  velocity 
perturbation  is  subtracted  from  the  lateral  velocity  perturbations.  After  about  five  meters  the 
slope  of  the  curves  becomes  essentially  constant.  This  constant  decrease  in  this  region  can  be 
attributed  to  the  viscous  effect. 

Calculated  entropy  production  rates  deviated  significantly  from  the  expected  values. 
As  mentioned  in  the  results,  the  viscosity  value  used  in  equation  6.2  may  not  be  accurate. 
Further  investigation  using  the  effective  viscosity  in  the  turbulence  cases  will  be  performed. 
At  present  AVUS  does  not  have  capability  for  this  calculation.  However  from  other  studies 
concerning  2D  airfoils  it  has  been  seen  that  the  effective  viscosity  can  be  as  much  as  100  to 
400  times  higher  than  the  natural  viscosity  in  some  cells.  This  would  place  our  entropy 
predictions  close  to  similar  order  of  magnitude  as  the  expected  values. 

A  computational  domain  in  any  CFD  simulation  cannot  be  infinite.  When  modeling 
airfoils  and  finite  wings  it  is  desirable  to  have  the  top,  bottom,  front,  and  side  faces  in  the 
freestream  region.  Having  the  back  face  in  the  freestream  region  may  require  too  many  grid 
points.  Therefore  it  may  be  possible  that  the  wake  and  trailing  vortex  sheet  may  extend 
beyond  the  back  face.  Thus  the  volume  integration  may  not  capture  all  of  the  entropy 
production  caused  by  the  wing.  Depending  on  how  much  entropy  is  being  produced  in  the 
wake,  this  could  cause  the  entropy  production  values  calculated  from  the  CFD  simulations  to 
deviate  significantly  from  the  expected  values.  When  comparing  two  wings,  however,  this 
should  not  matter  as  long  as  the  domains  are  of  equal  shape  and  volume. 
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CONCLUSION 


We  conducted  this  study  over  the  summer  of  2004.  A  technique  to  create  twisted 
wing  geometries  to  give  a  desired  lift  distribution  was  demonstrated.  A  grid  to  capture  the 
important  flow  physics  was  generated  on  the  twisted  wings.  With  this  gridding  methodology 
we  obtained  the  accurate  lift  values,  but  and  the  drag  predictions  were  shown  to  be  only 
reasonable.  We  have  also  demonstrated  the  method  of  Trefftz  plane  to  calculate  the  induced 
drag.  This  method  resulted  in  fairly  accurate  results  for  induced  drag. 

Some  interesting  results  were  obtained  from  the  CFD  simulations.  One  of  them  is 
that  the  parabolic  case  gave  less  induced  drag  that  than  the  elliptic  case.  This  goes  against 
the  classical  lifting  line  theory  which  says  the  elliptic  lift  distribution  yields  the  minimum 
induced  drag.  Also  the  parabolic  case  gave  a  lower  entropy  production  than  the  elliptic  case. 
This  agrees  with  Greene’s  [9]  hypothesis  that  says  that  parabolic  lift  distribution  will  yield 
minimum  entropy  generation.  However,  it  is  important  to  note  that  the  solutions  are  neither 
fully  converged  nor  grid  independent.  The  actual  lift  distributions  were  not  verified  and  will 
be  a  focus  of  future  work.  At  best,  we  have  developed  the  necessary  methodology  to  proceed 
into  the  work  to  follow. 
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FUTURE  WORK 


Since  computational  restraints  hindered  grid  quality  and  cell  count,  grid  independence 
studies  were  not  possible.  This  will  be  the  focus  of  future  work.  We  will  conduct  further 
CFD  simulations  using  Fluent  which  allows  for  solution-adaptive  grid  refinement  by  adding 
cells  in  the  area  of  interest. 

Fluent  also  has  the  ability  to  allow  the  user  to  add  cells  anywhere  in  the  domain.  This 
feature  should  help  improve  grid  quality  in  areas  that  have  large  jumps  in  cell  size. 
Performing  this  task  in  Fluent  allows  the  user  to  add  grid  points  in  a  specified  region  and  not 
affect  any  other  region  within  the  domain.  This  allows  for  efficient  grid  refinement  and 
should  help  with  solution  convergence  which  was  not  achieved  during  the  summer  of  2004. 

We  will  also  attempt  to  evaluate  the  actual  lift  distribution  that  the  wings  are 
producing  to  see  the  validity  of  the  lifting-line  theory  and  our  model  generation  technique. 

To  give  more  than  two  data  points  linking  drag  with  exergy  destruction,  a  wing  with 
an  arbitrary  lift  distribution  will  be  modeled.  This  wing  will  be  constructed  to  give  the  same 
total  lift  and  will  be  composed  of  the  same  airfoil  sections  used  in  this  study.  Doing  this  will 
provide  more  insight  into  the  connection  between  drag  and  exergy  destruction. 

We  expect  that  the  future  work  will  yield  grid  independent,  converged  solutions  of 
the  cases  performed  thus  far.  From  this,  the  drag,  entropy  production  and  exergy  destruction 
rates  will  be  related  in  terms  of  thermodynamic  quantity.  Using  these  results,  we  can  then 
proceed  to  validate  the  classical  lifting-line  theory. 
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Figure  6.4(a):  Rolled  up  streamlines  for  elliptic  case 


Figure  6.4(b):  Rolled  up  streamlines  for  parabolic  case 
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Figure  6.5:  Total  entropy  generation  rate  for  elliptic  turbulent  case  atz  =  0.5m 
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Figure  6.6:  Total  entropy  generation  rate  for  elliptic  turbulent  case  at  z  =  1.5m 
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Figure  6.7:  Total  entropy  generation  rate  for  elliptic  turbulent  case  at  z  =  2.9  m 
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Figure  6.8:  Total  entropy  generation  rate  at  leading  edge  for  elliptic  turbulent  case  atz  =  2.9 
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Figure  6.9:  Total  entropy  generation  rate  at  trailing  edge  for  elliptic  turbulent  case  atz  =  2.9 
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Figure  6.10:  Total  entropy  generation  rate  for  elliptic  case  at y  —  0m  (top  view). 
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Figure  6.11:  Total  entropy  generation  rate  for  parabolic  case  at  z  =  0.5m 
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Figure  6.12:  Total  entropy  generation  rate  for  parabolic  case  at  z  =  1.5m 
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Figure  6.13:  Total  entropy  generation  rate  for  parabolic  case  at  z  =  2.9m 
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Figure  6.14:  Total  entropy  generation  rate  for  elliptic  case  at  y  =  Om 


Figure  7.1:  Induced  drag  coefficients  at  varying  streamwise 
locations  with  vnrvim*  Trefftz  nlane  area  ( l;l Untie  Turbulent  Case! 
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cDi 

Cd  (Surface  Integration) 

Lifting-Line 

Trefftz  Plane 

Elliptic  Laminar  NS 

0.009839 

0.007909 

0.027967 

Elliptic  kco 

0.008127 

0.015342 

Parabolic  Laminar  NS 

0.010956 

0.006933 

0.016550 

Parabolic  kco 

0.007754 

0.017092 

Table  6.1:  Induced 

[  and  total  drag  coefficients 

CL 

Lifting-Line  Theory 

Surface  Integration 

Elliptic  Laminar  NS 

0.43064 

0.42644 

Elliptic  kco 

0.40499 

Parabolic  Laminar  NS 

0.43064 

0.44475 

Parabolic  kco 

0.43673 

Table  6.2:  Induced  and  total  drag  coefficients 


Elliptic 

Parabolic 

Laminar 

Turbulent 

Laminar 

Turbulent 

Cd, press 

2.7413xl0~2 

1.11 52x1 02 

1.2602x1  O'2 

1.262  lxl  O'2 

Cd,  friction 

5. 5440x1  O'4 

4.1898xl0‘3 

3.9471xl0'3 

4.4766xl0'3 

Cd, total 

2.7967xl0'2 

1.5342xl0"2 

1.6550xl0'2 

1.7092x1  O'2 

Cl, press 

4.2627xl04 

4.0516X101 

4.4489xl0_1 

4.3689xl0_1 

Cl,  friction 

1.7429xl0‘4 

-1.71 16xl0'4 

-1.3489x1  O'4 

-1.6052xl0"4 

Cl, total 

4.2644x1 0"1 

4.0499X10'1 

4.4475xl0_1 

4.3673X10'1 

Table  6.3:  Drag  and  lift  coefficients  calculated  with  surface  integration  over  the  wing  surface 
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Parabolic 

Elliptic 

Total  Sgen  Rate 

Laminar 

Turbulent 

Laminar 

Turbulent 

Expected  Values  (W/K) 

66.4 

68.6 

112.2 

61.56 

CFD  Results  (W/K) 

0.63758 

0.42924 

1.23474 

0.46412 

Table  6.4:  Expected  total  entropy  generation  rate  and  the  results  from  the  CFD  simulations 
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